EFFICIENT SUPERIMPOSITION RECOVERING ALGORITHM 



Han Li* Kun Gai*^ Pinghua Gong* Changshui Zhang* 

* Dept. of Automation, Tsinghua University, Beijing, China 
t Alibaba Group, Beijing, China 



t— i ■ 
O 

<n: 

>: 
o. 

£: 

OV 



> 

U 

o 



> 

o 
m 
^1- 



X 



ABSTRACT 

In this article, we address the issue of recovering latent trans- 
parent layers from superimposition images. Here, we assume 
we have the estimated transformations and extracted gradients 
of latent layers. To rapidly recover high-quality image lay- 
ers, we propose an Efficient Superimposition Recovering Al- 
gorithm (ESRA) by extending the framework of accelerated 
gradient method. In addition, a key building block (in each 
iteration) in our proposed method is the proximal operator 
calculating. Here we propose to employ a dual approach and 
present our Parallel Algorithm with Constrained Total Vari- 
ation (PACTV) method. Our recovering method not only re- 
constructs high-quality layers without color-bias problem, but 
also theoretically guarantees good convergence performance. 

Index Terms — Superimposition recovering, proximal 
operator, optimization 



Assume we have the extracted gradients and transforma- 
tions of each latent transparent layer, the reconstruction step 
is the final crucial part for reflection separation. We assume 
the variations of transmitted layers in each mixtures conform 
to a parametric transformation f(x,9) (x is the pixel coor- 
dinates) with different parameters 0j. Here we propose an 
Efficient Superimposition Recovering Algorithm (ESRA) to 
fast recover the high quality latent layers. 

1. EFFICIENT SUPERIMPOSITION RECOVERING 
ALGORITHM 

With estimated transformation parameters 9i, we align the 
transmitted layers by warping mixtures Ii with f^ 1 (x,9i). 
Then our mixing model is rewritten as: 

ktfi 1 ®) = aaLtty+aaL^ifr 1 ®), i = 1,... , m . 

(1) 

Here L* is the latent transmitted layer, L r ® is the reflected 
layer in jfh (mixtures), an,ai2 is the mixing coefficients. 
With this new mixing model, the influence of parametric 
transformations f(x,9i) can be ignored in the interme- 
diate recovering process. For simplicity, we use Ii(x) to 
represent Ii(f~ 1 (xj). Li(x) and Li + i(x) denote and 



a2iL r ®(f~ 1 (x)), respectively. Let E l (x) stand for the ex- 
tracted gradients from L l (x). To recover high quality latent 
image layers, we propose to employ L\ penalty on the ex- 
tracted gradients and nonnegative constraints on the layers' 
intensities along with the L2 loss of the mixing model. Thus 
our recovering objective function is written as: 



m+l 

min F(l vec ) = X y \WL t (x)-E t (x) 

0<I'«<1 ^ 1 



(2) 
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'a i iLi(x)~L j+ i(x)f 



where l l 



{L 1 ), ■ ■ ■ , vec T (L m+1 )] T is a large vec- 
tor containing all pixel values in all latent layers. The first L\ 
term enforces the agreement between reconstructed layer gra- 
dients and extracted layer gradients, while the second L2 term 
tends to satisfy our mixing mode. Since the extracted gradi- 
ents are nonzero at very few coordinates, the L\ norm term 
not only prefers layers with sparse gradients but also avoids 
over-smooth results. A is a trade off coefficient. 

To solve the nonsmooth convex optimization model (O 
efficiently, we denote 

^(i^-anL^-L^f, s.tO < l vec < 1, 



f(l l 



E 



m+l 



5 (r ec )=A ^ \VL\x)-E\x)\. 

X, 2 — 1 

(3) 

Here g(l vec ) is the £1 penalty on the extracted gradients and 
f(l vec ) corresponds to the L2 loss and nonnegative con- 
straints. f(l vec ) can be formulated in the following matrix 
form: 

1, 



f(i vec ) = ~\\Ar 



b\\ 2 , s.t. < r ec < 1, 









"vec(Ji) 


where A = 




,b = 






_a m il I 




_vec(/ m ) 



(4) 



where f(l vec ) is continuously differentiable and V f(l vec ) — 
A T (Al vec ~b), of which Lipschitz constant L(f) = \ max {A T A) - 
J2t afi + !- and 1 € R hwxhw is the unit matrix. We note the 
objective function in (fJJ is a composite function of a differen- 
tial term f(l vec ) and a non-differential term g(l vec ). Denote 



m+l 



■l|r 



(5) 



which is the first order Taylor expansion of f(l vec ) at l^x, 
with the squared Euclidean distance between l vec and Vj^i as 
the regularization term. The traditional gradient descent al- 
gorithm obtains the solution at the fc-th iteration (fc > 1) by 
l v k ec = argminPi si i«=c i (Z , ' ec ) + g(l vec ) with a proper step 
size L s (greater than £■(/)). Here we propose to employ the 
accelerated gradient descent (HE] to solve the reconstruction 
problem, named Efficient Superimposition Recovering Algo- 
rithm (ESRA). Here we generate a solution at the fc-th itera- 
tion (fc > 1) by computing the following proximal operator 



JVC 

'ft 



a rg Q T^ <i p L „Y k (r" sc ) + g(i uec ) (6) 



where Y x = l v ec and Y k = l v k e _\ + **-^~ x - Z£I C 2 ) for 
fc > 1. We note that Yk is a linear combination of lY^ and 
The combination coefficient plays an important role in 
the convergence of the algorithm. As suggested by 0, we set 

t = 1 and t k = (1 + -y/ifc_i + l)/2 for fc > 1. Accord- 
ing to the theoretical analysis in [3], this accelerated gradient 
descent method can get within 0(l/fc 2 ) of the optimal ob- 
jective value after fc steps. While solving problem © is still 
very challenging, we propose a Parallel Algorithm with Con- 
strained Total Variation (PACTV) method to find the optimal 
solution, which is presented in the sequel. 

Algorithm 1: ESRA(L S , N) 

Input : L s > L(f)- An upper bound on the Lipschitz 
constant of V/. N is the total number of 
iterations. 

1 StepO: Take Y t = Z uec = (Mro+1)il)l t a = 1. 

2 Stepk:(fc > 1) Compute 

ll ec = P LttYh {l^)+g{l^), (7) 



t k ■■ 

z 

\r ivec | tk— 1 ~ 1 

Y k+1 - tfc H t 

Eft 



(/r c - ra 



(8) 
(9) 



Output: Z)(f c is the final recovered result. 



ll ec = argmin { £ £ (A|V V (x) - & (a 



i—1 x 



0<J£ ec <l 

+ ^|| L i( a .)_ di ( a .)|P)}. 



(10) 



As illustrated in ( [Tol l, finding Z£i ec is to solve following m + l 
separable problems with constrained total variation in paral- 
lel: 

of&E {l\\L(x)-d(x)\\ 2 + P\VL(x)-E(x)\). 

Here /? = X/L s , and L,d,E represent L l ,di,E l , respec- 
tively. Similar with the image denoising problem HO, we 
propose a dual approach to solve ( fTTT i and give some notation 
in order: 

• V is the set of matrix-pairs (p, q) where p € R(' 1_1 ) X1 " 
and q € ]R to < ,i, - L ) that satisfy \p id \ < 1 and \q id \ < 
1 And we assume p ,j = Ph,j = Qi,o = Qi,v> = 
0, for every i = 1, • • • , h, j = 1, • • • , w. 

• The linear operation C : R^- 1 )™ x K^^i) -> R hxw 
is defined by the formula C(p, q)ij = Pi-ij + qi,j-i — 
Pi,j+Qi,j Vi,j. 

• The operator C T : R hxw -> E^ 1 ^ x RM«^-i) which 
is adjoint to £ is given by C T (L) = (p,q), where 
Pi,j ~ ^i+ij' — and qij = — ^i.j- 

• Pc is the orthogonal projection operator on the convex 
closed set C = {L : < L < 1}. 

Equipped with these notation, we derive a dual problem of 
( fTTT i. and give following proposition to state the relation be- 
tween the primal and dual optimal solutions. 

Proposition 1. Let (p, q) € V be the optimal solution of the 
problem 

min {H(p,q) = h-\\H c {d- (iC{p,q))\\ 2 + 

2 (12) 

\\d-f3£(p,q)\\ 2 ) +(3[Tr(p T E 1 ) + Tr(q T E 2 )] }. 

where H C {L) = L - P C {L) for every L E R hxw . Then the 
optimal solution of Uli is given by L — Pc{d — /3£(p, q)). 



2. PACTV VIA DUAL APPROACH 



Proof First note the following relation holds true: 



Given problem (O, we observe it can be solved block sepa- 
rable in the following way. If we denote Yk — j^-S7f(Yk) = 



vec 



(di),- 



, vec 



(dm+i)} (di e 



vhxw 



we can split Yk— -r-Vf(Yk) into ; 



i = 1,- • ■, m + l), 
1 separable parts. Then 



by employing the definition of y), we transform (O into the 
following form: 



max{px : |p| < 1}. 
p 



Hence, we can give 



ft 



\V k L-E k \ = max T(L,p,q), 
( P ,q)ev 



(13) 



(14) 



where, 



p hxw 



defined by 



h-l w-1 

T(L,p,q) =J2Y1 [Pi,j( L md - L i,j ~ E uJ 
i=i j=i 

+ - - E 2i J )] 



h-i 



(15) 



+ 7,Phj( L h,i±l - Lh,j - E 2h j ). 
3=1 

With this notation we have 

T(L,p, q) = Tr(£(p, q) T L) - Tr(p T E,) - Tr(q T E 2 ) 

Thus the original problem ( fTTT i becomes 

i max . {hL-d\\ 2 + p[Tr(£(p,q) T L) 



(16) 



mm 

0<L<1 (p, q )j 



Tr(p T E x ) - Tr(q r E 2 )] }. 



(17) 



Since the objective function is convex in L and concave in 
p. q, we can exchange the order of the minimum and maxi- 
mum and get 



(18) 



max mm 

{p,q)ev 0<L< 



max min - d\\ 2 + /3\Tr(£(p, q) T L) 

(p,q)£V0<L<l 12 L \ \ J J 

— Tr{p T E\) — Tr(q T E 2 )~\ |. 
and which can be written as 

i^ {1 [\\L - (d - (3£(p, q))\\ 2 - \\d - /3C(p, q)f 

+ \\d\\ 2 ]-(3[Tr(p T E 1 )+Tr(q T E 2 )]}. 

(19) 

Thus the optimal solution of the inner minimization problem 
is 

L = P {0 <L<i}{d-PC(p,q)). (20) 

And last, we plug the above expression for L back into ( fT9l 
and ignore the constant term, we obtain the dual problem is 

min {H{p,q) = \(-\\Hc(d-(3C{p,q))\\ 2 + 

\\d - 0£(p, q)\\ 2 ) + /3[Tr(p T E 1 ) + Tr(q T E 2 )}}, 

which is the same as ( fT2l . □ 

what's more, given (fT2l . we can easily have following 
lemma. 

Lemma 1. The objective funtion H of ( 1721) is continuously 
differentiable and its gradient is given by 
VH(p,q) = -(3£ T P c (d~/3£{p,q)) + p(E u E 2 ). (21) 

And let L(H) be the Lipschitz constant of V H(p, q), then 
L(H) < 8/3 2 . 



Proof. Consider the function s 

s(L) = \\H C {L)\\ 2 . (22) 
Then the dual function (fT2l can be written as: 

H(p, q) = \{-s{d- pC{ Pl q)) + \\d - pC(p,q)\\ 2 ) 
+(3[Tr(p T E 1 )+Tr{q T E 2 )]. 

Obviously, s(-) is continuously differentiable and its gradient 
is given by 

Vs(L) = 2(L - P C {L)). (24) 

Therefore, 



s(d - PC(p,q)) + \\d- pC{p, q)\\ 2 ) + (3(E U E 2 ) 



-In- 

= ±/3£ T (Va(d - PC{ P , q)) - 2(d - 0C(p, <?))) + /3(E 1 ,E 2 ) 
= -l3£ T P c (d - (3£(p, q)) + p(Et, E 2 ) 

(25) 

Then for every two pairs of matrices (pi, q±), (p 2 , q 2 ) 
where Pl € and q t e l^f^ 1 ) for i = 1, 2, we have 

||Vfl-(pi,gi)-Vfl-(pa,«&)|| 

= fi 1 1 £ T [P C (d - pCfa , Ql ))]- £ T [P C (d - (3£{p 2 , q 2 ))] 1 1 
< f3\\£ T \\\\P c (d - f3£( Pl , qi )) - P c (d - p£( P2 ,q 2 ))\\ 
<(3 2 \\£ T \\\\£(p 1 ,q 1 )-£(p 2 ,q 2 )\\ 

<p 2 \\£ T \\\\£\\\\( Pl ,q 1 )-(p 2 ,q 2 )\\ 

= /3 2 \\£ T \\ 2 \\(p 1 ,q 1 )-(p 2 ,q 2 )\\ 

(26) 

here the above inequalities follow from the non-expensiveness 
property of the orthogonal projection operator and prop- 
erty of linear operators £,£ T . And from [4], we have 
\\£ T (x)\\ < V8\\x\\. Therefore, implying that ||£ T || < \/8 
and hence L{H) < 8/3 2 . □ 

With definition of H(p,q) and WH(p,q), fast gradient 
projection (FGP) is applied on the dual problem (fT2l . And 
the complexity of each iteration in FGP is O(hw). Above 
all, our proposed Parallel Algorithm with Constrained Total 
Variation (PACTV) is using FGP to solve the rn + 1 dual 
problems ( fT2l in parallel. Then we catenate the optimal L* 
(i = 1, . . . , m+1) and resize them into vector form to achieve 

jvec 
L k ■ 

Given above proposition and lemma, we can use the fast 
gradient projection (FGP) on dual problem ( fT2l . Fast gradient 
projection (FGP) is outlined in Algorithm [2] Here P-p(p,q) 
means projecting the matrix-pair (p, q) on the set V . And 
finally we achieve the optimal solution of dTTb . Then our re- 
covering method ESRA is outlined in Algorithm^ 

In our implementations, we set the total iteration number 
of ESRA is 100 and FGP tolerance is 0.0001, and we also 



Algorithm 2: FGP(6, /3, N, E x , E 2 ) 



Input : d e R hxw , = X/L s , N is the total number 
of iterations, E%,E2, 

1 Step 0: Take 

ipuQi) = (po,Qo) = (0(^i)x«;,0hx(«^i)))*i = !• 

2 Step k:(n > 1) Compute 

{Pk,qk) = Bp[(Pk,qk) - ^V#(pfc,g fe )], (27) 

1 + ^/l + 4t£ 
ifc+i — ^ , (28) 

(pfcfij^fefi) = (Pfc)9fc) + (-7 — -)(Pk-Pk-x,qk-qk-i) (29) 

Output: L* An optimal solution of ( fTTT i up to a 
tolerance. 



set L s — 2L(f) to ensure a constant stepsize. The initial 
value of l vec is zero. The final recovered reflected layers of 
(f2|i should be warped with fa and enhance the intensity by 
2 to be visible. Our recovering method launches a general 
optimization framework and can be extended to solve other 
reconstruction problems in J5] |6) . 
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